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Abstract 



In an attempt to explain the uniqueness of the coding mechanism of Hving 
cells as contrasted with the multi-species structure of ecosystems we examine 
two models of individuals with some replicative properties. In the first model 
the system generically remains in a multi-species state. Even though for some 
of these species the replicative probability is very high, they are unable to 
invade the system. In the second model, in which the death rate depends on 
the type of the species, the system relatively quickly reaches a single-species 
state and fluctuations might at most bring it to yet another single-species 
state. 
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I. INTRODUCTION 



The problem of the origin of hfe and its early evolution is clearly one of the most fun- 
damental problems of modern science. In spite of considerable efforts, even most basic 
questions in this multidisciplinary issue remain unanswered. The basic frame of most the- 
ories of the emergence of life was set many years ago by Oparin [|l|, who proposed that 
life emerged as a result of gradual evolution of nonorganic matter. Oparin's ideas were to 
some extent verified by experiments done by Miller and Urey [Q], who showed that for a 
system of water and an atmosphere consisting of gases that were thought to be common 
on the pre-biological Earth, electrical discharges resulted in the formation of some amino 
acids and nucleotides. It is believed that once created in sufficient concentrations, these 
molecules entered complicated synthetic reactions, which produced more and more complex 
molecules. Some of these complex molecules had catalytic and presumably to some extent 
even autocatalytic properties The autocatalytic molecules (or rather systems of them), 
if of sufficient stability, were clearly more likely to survive in a competition for reactants. 
The gradual evolution of such autocatalytic systems, subjected to Darwinian selection, re- 
sulted in mastering their surviving skills and eventually led to the emergence of life. Even 
though the above-sketched scenario might seem plausible, many of its important details are 
still unresolved |Q. 

The replicatory mechanism of contemporary living cells is very sophisticated and basi- 
cally unchanged since the emergence of the first living cells, which presumably took place 
about 3.5-4 billion years ago. One of the characteristic features of this mechanism is that 
the tasks of coding and catalysis are being assigned to different macromolecules, namely to 
nucleic acids and proteins, respectively. Moreover, the code, i.e., the way amino acids are 
encoded by nucleotides, is universal for all living cells. It suggests an interesting possibility 
that all living cells are actually descendants of a single pra-cell, which happened to develop 
the most effective surviving skills. 

At the same time, however, it raises some questions. One might expect that in search 
for the most effective cells Nature tried many variants of different effectiveness. Why did a 
certain code predominated all other variants? Was this variant really of such an enormous 
effectiveness or maybe predomination was somehow a generic feature of prebiotic dynamics? 
One can notice, however, that in any contemporary ecosystem a large number of species 
coexist and these species are clearly of different effectiveness. Nevertheless the invasion of an 
ecosystem by a single species is an unobservable phenomenon. Although very much different, 
both contemporary ecosystems and primeval soup might be regarded as composed of certain 
individuals with some replicative properties. Why thus did Nature select the single- species 
solution at the early stages of life and why does it prefer multi-species solutions at later 
stages? 

Various aspects of the problem of the origin of life and biological evolution have been 
already modeled 0. Even though such models are, by necessity, highly simplified, they help 
us to understand the essence of these complex phenomena. For example, one can construct 
simple models of biological evolution which explain why the dynamics of extinction of species 
has some scale-invariant properties [Q. 

The problem of multiplicity of species in replicative systems has been also already ad- 
dressed in the literature. A comprehensive review of bio-chemical aspects of this problem 
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was written by Orgel [0 and most recently by Szathmary [§. Certain simple models of 
replicative systems have been also recently examined For example, in the model dis- 

cussed by Szathmary and Maynard Smith an ensemble of replicators is described in terms 
of differential equations. In particular, the concentration of the i-th replicator is described 
by the following equation p]: 

dxi/dt = kiX^, (1) 

where ki is the growth rate constant of the i-th replicator and p is the order of replication. 
It turns out that asymptotic (t oo) concentrations (or rather their ratios) depend on p. 
For p = 1 (Malthusian growth), the replicator with the largest ki becomes dominant and 
such a case is characterized as 'survival of the fittest'. For p < 1 (parabolic growth), the 
ratios of concentrations become finite, which is termed as 'survival of everybody'. In the 
case p > 1 (hyperbolic growth), the dominant replicator is the one with the largest product 
of the initial concentration and the growth rate constant, which is termed as 'survival of 
the common'. The interest in growth laws with p 7^ 1 is partially motivated experimentally, 
since it was shown that certain oligonucleotides, which presumably played an important 
role in prebiotic dynamics, indeed follow the growth law with p < 1 |]TT],|T2[. Since such a 
growth implies 'survival of everybody', we are faced with a problem of the transition to the 
Malthusian growth, which would explain 'survival of the fittest'. However, recently it was 
shown by Lifson and Lifson that for more general growth laws than Eq. (|l]), the 'survival of 



the fittest' takes place even in the p < 1 case |jTO| 



An important assumption underlying models leading to differential equations like Eq. (|I]) 
is that replicators are perfect, i.e., a replicator produces at a certain speed its exact copy. 
In our opinion, to model evolution at early stages we should rather consider a system of 
imperfect replicators, which, for example, would produce their copies only with a certain 
probability and otherwise they would mutate. 

In the present Paper, we examine two simple models of such systems. In our models. 



replicators, which might replicate or mutate, exist in infinitely many varieties |T^. As our 
main result, we show that behaviour of these models strongly depends on some details of 
dynamics of these models. Namely, only in one of these models the evolution proceeds in a 
single-species way, i.e., with the majority of the system descending from the same ancestor. 
In the second model, such single-species states are very unstable and the system evolves 
through multi-species paths. 

In section II we define our models and present their basic properties. In section III we 
examine in more details the behaviour of each model, emphasizing the difference between 
single- and multi-species evolution. Section IV contains our conclusions. 

II. MODELS AND THEIR BASIC PROPERTIES 

A. Model I 

Let us consider a system composed initially of L individuals. These individuals might 
be regarded as complex molecules at the prebiotic era immersed in the primeval soup and 
thus involved in a number of catalytic or autocatalytic reactions. With each individual we 



3 



assign randomly the replication probability pi {0 < pi < 1 for i = 1,2 ... ,L) that the i-th 
individual will exactly replicate itself in the course of reproduction. The dynamics of this 
model, which in the following will be referred to as Model I, is specified as follows: 

1. Choose an individual at random. The chosen individual is denoted by i. 

2. With the probability L/N the i-th individual dies. The constant ^ 1 might be 
regarded as a certain 'environmental capacity'. Namely, provided that initially we 
have L < N, L will never exceed A^. 

3. With the probability 1— L/N the i-th individual survives and produces a new individual 
j. The probability pj assigned to the j-th individual is equal to pi (parent's value) 
with the probability Pi and is equal to a random number from the interval (0,1) with 
the probability 1—pi. This rule means that if a copying error happens, it dramatically 
changes the properties of the new individual. 

These rules imply that in the steady state the death rate {L/N) equals to the reproduction 
rate (1 — L/N) and thus on average L = N/2. As far as the number of individuals is 
concerned. Model I is equivalent to a certain random walk problem. Indeed, from the rules 
stated above one can infer the following equation for the probability P{L,t) that in our 
system there are L individuals at time t: 

P(L, t + A) = ^P(L + 1, t) + (1 - ^)P(L - 1, t). (2) 

This equation describes changes in our system after updating a single molecule. To conform 
to the Monte Carlo simulations presented later, we assume that such a single update takes 
J of time (i.e., a unit of time corresponds to a single, on average, update of each individual). 
Thus, in Eq.(|^) we have A = Introducing the variable M = L — we can rewrite Eq.(^ 

as 

A^ 1 M N 1 M N 

P{- + M, t + A) = (- - -)P(- + M - 1, t) + (- + -)P(- + M + 1, t), (3) 

which is clearly the equation of a random walk with attraction toward M = 0, i.e., L = 
Y- It has been already shown that the so-called 'dog-flea' model is also equivalent to 
a similar random walk problem with attraction and that fluctuations in this model 
around equilibrium (M = 0) become negligible in the limit, which in our case corresponds 
to A^ oo. Thus, we expect that in Model I fluctuations of the number of individuals 
around L = ^ are also small for large A^. 

To examine replicative properties of our model, we resort first to Monte Carlo simulations. 
Since the implementation of the above rules on the computer is rather straightforward, we 
present only the results of these simulations. In Fig.|l] we present the time evolution of 
the average replication probability p{t) = j^Y.i=iPi- Simulations were done for A^ = 10^ 
and initially the probabilities pi {i = 1, . . . , N) were chosen at random. In all simulations 
reported in the present paper the initial number of individuals L is equal to A^. One can 
clearly see that the average replication probability increases in time but the increase is very 
slow and it is not obvious what value is reached in the steady state (i.e., for t — > oo). 



4 



However, below we present some analytical calculations which show that if the limit 
— > cx) is taken first, then for t — oo the average replication probability p{t) converges 
to unity. Our strategy is to write the evolution equation for the higher order moments 
of replication probability and then to solve the resulting infinite set of equations in the 
steady state. First let us assume that the evolution in our model lasted long enough to 
equilibrate it with respect to the number of molecules L. Thus, we approximate the death 
and reproduction probabilities as^ = l — ^ = |. Introducing the notation p'(t) = j J2i=iPi, 
where / = 1, 2 . . . , we can write the following evolution equation for the first moment of Pi 
(i.e., for p{t)) 

jP\t + A) - ^p\t) = -ipi(t) + l\p\t) + 1(1 - p\t))], (4) 

were we put y as an average number of individuals. The first and second terms on the 
right-hand side of Eq.(^) describe the changes due to a single update caused, respectively, 
by the death and reproduction processes. The term |(1 —p^{t)) corresponds to production 
of an individual with a randomly assigned replication probability, which thus on average 
takes the value ^. We can write similar equations for arbitrary I: 

jp\t + A)~^p\t) = -lp\t) + l[p^^\t) + j^{l-p\t))], for/ = l,2... . (5) 

In Eq.(^ we used the fact that the l-th moment of a random variable, which is uniformly 
distributed on (0,1), is equal to Jq s''ds = j^. In the steady state, the left-hand side of 
Eq.(^ is zero and thus we obtain 

=p'+i + -1-(1 -pi) for Z = l,2... . (6) 

I ~\~ -L 

This infinite set of equations can be solved. Namely, when we add the Eqs.(^) for Z = 1, 2, . . . , 
all higher-order moments cancel out and we obtain 

OO -| 

p' = a-i)E-^ (7) 

i=2 ^ 

and thus p^ = J2'^2 \l I^i^i j — ^i since both series diverge and differ only by unity. Sim- 
ilarly, all other moments in the steady state are equal to one. Thus we expect that in 
the limit t ^ oo the replicative probability p{t) in Fig.|I] increases to unity even though the 
convergence seems to be very slow. 

B. Model II 

Before discussing other properties of Model I, let us consider the model where the prob- 
ability of death of a certain individual depends not only on the total number of individuals 
(as in Model I) but also on the individual itself. Such a modification is motivated by the fact 
that in the primeval soup the survival of a molecule was determined not only by the access 
to substrates (and then the total number of molecules is likely to determine the death rate) 
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but also by the stability of a given molecule against, e.g., radiation, which clearly depends 
on the type of this molecule. 

Thus, let us consider the model which in the following will be referred to as Model II. 
To each individual i, in addition to the replication probability Pi we assign randomly certain 
individual survival probability rj (0 < < 1). The dynamics of this model is specified as 
follows: 

1. Choose an individual at random. The chosen individual is denoted by i. 

2. With the probability L/N the i-th individual dies due to the lack of reproductive 
substrates. 

3. Provided that the individual survived the previous step: 

• it dies with the probability 1 — rj, 

• it survives with the probability and reproduces according to the rule analo- 
gous to that of Model I. Namely, the new individual with the probability pi has 
the same replication probability and death probability as its parent (i.e., Pi and 
Tj, respectively) and with the probability 1 — Pi these probabilities are chosen 
randomly anew. 

For Ti = 1 {i = 1, . . . , L) Model II becomes equivalent to Model I. Monte Carlo simu- 
lations for Model II are also straightforward and we present only the results. In Fig.^ we 
present the average replication probability p{t) defined in the same way as for Model I. One 
can see that the convergence to unity is in this model much faster than in Model I. 

III. SINGLE- VERSUS MULTI-SPECIES EVOLUTION 

But there are more important differences between these models than the rate of conver- 
gence. Certain indication of different behaviour is seen in Fig.|^, where we plot 1 — p{t) as a 
function of time in a logarithmic scale. One can see that the late-time evolution of Model II 
proceeds in steps between which the system basically remains at the same level of p{t) and 
none indication of such a stepwise behaviour is seen for Model I. This stepwise behaviour 
suggests that Model II remains mostly in a single-species state with majority of individuals 
belonging to the same species. Individuals i and j belong to the same species if pi = pj and 
Tj = rj] for Model I the second condition is always satisfied. 

To confirm such a scenario we present in Fig.|^ snapshot configurations for Models I and 
II after 5000 Monte Carlo steps and = 10^. Indeed, after this simulation time Model II 
was brought to a single-species state with only few individuals of pi ^ 0.9996 . . . (they are 
not shown in Fig.^ since they all have pi < 0.99). On the other hand. Model I still remains 
in the multi-species state. 

We would like to point out that, of course, for finite there exists a finite probability 
that Model I can be brought into a single-species state (and for small N one can indeed see 
such a behaviour in simulations) but for large (~ 10^) this would require an extremely 
long simulation time. What is, however, more important is that for Model I single-species 
states are relatively unstable and this model mostly remains in multi-species state. Such a 
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behaviour is clearly seen in Fig.^, where we show the time evolution of the percentage s of 
individuals belonging to a dominant species in the system. As an initial configuration we have 
chosen a single-species state with prescribed values of replication and death probabilities, 
namely we set Pi = ri = pq, where po = 0.9 or 0.99 for i = 1,2, . . . , N (for model I, rj = 1 
independently on po). One can see that Model I indeed quickly abandons the single-species 
state. We have checked that also for a larger po the behaviour of Model I is basically the 
same and the model quickly evolves toward the multi-species state. Evolution of Model II is 
different. When prepared in a state of large po (e.g., 0.99), the model remains in this state 
for very long time. Even when po is smaller (e.g., 0.9), this model, after some short transient, 
ends up in a single-species state. Further evolution of Model II consists of small fluctuations 
within such a state, which sometimes might be strong enough to bring the system, again via 
some short transient, to another single-species state and usually with larger replication and 
survival probabilities. 

Such a scenario is also confirmed in Figs.|^0. We simulated the system of the size 
N = 2 ■ 10^ and with random initial probabilities Pi and r^. In Fig.|^, which shows the 
occupation of dominant species s, one can see that Model I has rather irregular behaviour. 
Sometimes a dominant species occupies a great majority of the system {s ~ 1) but sometimes 
it is only a small fraction of the system. The behaviour of Model II is different (Fig.|^). The 
dominant species almost always occupies most of the system. Only during very rare and 
short-lived fluctuations, s becomes substantially smaller than unity. In Fig.^ we show the 
replication probability of the dominant species for the Monte Carlo runs shown in Fig.^ and 
Fig.p. One can see that in Model II, contrary to Model I, the dominant species are very 
longed-lived. A comparison of Fig.|^ with Fig.^ shows that even when the dominant species 
is unchanged, the percentage s occupied by it might fluctuate wildly. 

IV. CONCLUSIONS 

We have examined two very simple models of systems of replicative individuals. Al- 
though both models seem to evolve toward the state of perfect replicability, their evolution 
is markedly different. For Model II evolution proceeds through a sequence of consecutive 
transitions, between which the system remains basically in a single-species state (i.e, with 
almost all individuals being identical). In our opinion, this model might describe prebi- 
otic evolution until the invention of the universal code. According to this model, when the 
relatively stable and almost error-free replicative mechanism was found, it quickly invaded 
the whole system. On the other hand. Model I during its evolution remains mostly in a 
multi-species state. Such a multi-species structure resembles modern ecosystems, where a 
large number of species coexist and are constantly struggling for survival. 

Finally, let us notice that although the evolution of a system as a whole seems to proceed 
slower in Model I than in Model II (Fig.^, there are some species in Model I with replicative 
probabilities very close to unity (Fig.^ Fig.|^). Thus, multi-species dynamics in Model I 
enhances nucleation of species of very high effectiveness. However, even they are unlikely to 
'reign' for a long time. 
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FIGURES 



FIG. 1. The average replication probability p{t) as a function of time t for Model I (solid line) 
and Model II (dotted line). Simulations were performed for N = 10^ and initially the probabilities 
Pi, (i = 1, . . . , N) were chosen at random. 

FIG. 2. The logarithmic plot of 1 —p{t) as a function of t for N = 10^ and initial conditions as 
in Fig.|l[ 

FIG. 3. The replication probability pi as a function of i for Model I (small dots) and Model 
II (diamonds) after 5000 Monte Carlo steps and for A'^ = 10^. Almost all individuals for Model II 
have the same value of pi = 0.9996 . . . and the corresponding diamonds constitute a thick line in 
the upper part of the figure. 

FIG. 4. Occupation of a dominant species s as a function of time for Model I (bottom lines) 
and Model II (top lines). In the initial configuration all individuals are identical with the replication 
and death probabilities equal to 0.9 or 0.99. Simulations were done for N = 2 ■ 10^. 

FIG. 5. The occupation of the dominant species s as a function of time for Model I. Simulations, 
were done for N = 2 ■ 10^ and random initial probabilities pi (i = 1, 2, . . . , A^). 

FIG. 6. The occupation of the dominant species s as a function of time for Model II. Simulations 
were done for N = 2 ■ 10^ and random initial probabilities pi and rj (i = 1, 2, . . . , N). Three large 
fluctuations seen for t < 50000 resulted in changing the dominant species (see Fig.|^. 

FIG. 7. The replication probability pi of the dominant species for Model I (dots) and Model 
II (densely plotted o which constitute basically a thick line at pi ~ 0.997). 
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